#include "dyn_trapzDriver.h"

TrapzDriver::TrapzDriver(Functor &_func, double _a, double _b)
    : func(_func)
{
    setA(_a);
    setB(_b);
    n=0;
}

double TrapzDriver::computeNextRefinement()
{
    // We introduce some parameters.
    double x0;                      // Position of the first newly added interior point
    double totalNumPoints;          // Total number of interior points.
    double sum;                     // Sum of the contribution of the interior points.
    double delta;                   // Spacing between the new interior points.

    int interiorPoints,j;

    // We increment the refinement stage.
    n++;

    // We check to see if this is the first call.
    // If so, we compute the average.
    if (n==1)
    {
        // The s= assignment will always us to reuse that value in the future.
        return (s=0.5*(b-a)*(func(a)+func(b)));
    }
    else
    {
        // We compute the number of interior points.
        for (interiorPoints=1,j=1;j<n-1;j++)
            interiorPoints <<=1;    // Bitwise operation to double. Faster than doing *=2.

        totalNumPoints = (double) interiorPoints;   // Typecast for future use.
        delta = (b-a)/totalNumPoints;               // Spacing between the points.
        x0 = a+0.5*delta;
        for (sum=0.0,j=0;j<interiorPoints;j++)
        {
            sum += func(x0);         // We increment the sum.
            x0 += delta;             // And go to the next point.
        }

        // Something like a moving average.
        s = 0.5*(s+(b-a)*sum/totalNumPoints);
        return s;
    }
}
